function [channel_ir, channel_params, fig_impulse_response] = channel_impulse_response(fs, t_max, num_paths, channel_type)
% 信道冲激响应分析
% 输入参数：
%   fs - 采样频率 (Hz)
%   t_max - 最大时间 (秒)
%   num_paths - 多径数量
%   channel_type - 信道类型 ('urban', 'suburban', 'indoor', 'rural')
% 输出参数：
%   channel_ir - 信道冲激响应
%   channel_params - 信道参数
%   fig_impulse_response - 图形句柄

% 添加路径
addpath('../Common');
colors = color_definitions();

%% 并行计算设置
[use_parallel, pool_info] = parallel_manager();

% 时间向量
t = 0:1/fs:t_max-1/fs;
channel_ir = zeros(size(t));

% 根据信道类型设置参数
switch lower(channel_type)
    case 'urban'
        % 城市环境 (密集多径)
        max_delay = 2e-3; % 最大时延 2ms
        delay_spread = 0.5e-3; % 时延扩展
        relative_powers = linspace(0, -20, num_paths); % 功率衰减
        
    case 'suburban'
        % 郊区环境 (中等多径)
        max_delay = 1e-3; % 最大时延 1ms
        delay_spread = 0.2e-3; % 时延扩展
        relative_powers = linspace(0, -15, num_paths); % 功率衰减
        
    case 'indoor'
        % 室内环境 (较少多径)
        max_delay = 0.5e-3; % 最大时延 0.5ms
        delay_spread = 0.05e-3; % 时延扩展
        relative_powers = linspace(0, -12, num_paths); % 功率衰减
        
    case 'rural'
        % 农村环境 (很少多径)
        max_delay = 0.2e-3; % 最大时延 0.2ms
        delay_spread = 0.02e-3; % 时延扩展
        relative_powers = linspace(0, -8, num_paths); % 功率衰减
        
    otherwise
        % 默认城市环境
        max_delay = 2e-3;
        delay_spread = 0.5e-3;
        relative_powers = linspace(0, -20, num_paths);
end

% 生成随机时延
delays = sort(rand(1, num_paths) * max_delay);
delays(1) = 0; % 第一条径时延为0

% 生成信道冲激响应
if use_parallel
    % 使用并行计算生成多径分量
    fprintf('使用并行计算生成信道冲激响应...\n');
    
    % 预分配临时结果存储
    temp_results = cell(num_paths, 1);
    
    parfor i = 1:num_paths
        % 瑞利衰落系数
        alpha = (randn + 1i*randn)/sqrt(2) * 10^(relative_powers(i)/20);
        
        % 冲激响应
        delay_samples = round(delays(i)*fs);
        if delay_samples < length(t)
            temp_results{i} = struct('delay', delay_samples+1, 'alpha', alpha);
        else
            temp_results{i} = struct('delay', [], 'alpha', 0);
        end
    end
    
    % 将临时结果复制到最终数组
    for i = 1:num_paths
        if ~isempty(temp_results{i}.delay) && temp_results{i}.delay <= length(t)
            channel_ir(temp_results{i}.delay) = temp_results{i}.alpha;
        end
    end
else
    % 串行计算
    for i = 1:num_paths
        % 瑞利衰落系数
        alpha = (randn + 1i*randn)/sqrt(2) * 10^(relative_powers(i)/20);
        
        % 冲激响应
        delay_samples = round(delays(i)*fs);
        if delay_samples < length(t)
            channel_ir(delay_samples+1) = alpha;
        end
    end
end

% 信道参数
channel_params = struct();
channel_params.delays = delays;
channel_params.relative_powers = relative_powers;
channel_params.max_delay = max_delay;
channel_params.delay_spread = delay_spread;
channel_params.rms_delay_spread = sqrt(sum(delays.^2 .* 10.^(relative_powers/10)) / sum(10.^(relative_powers/10)));

%% 可视化结果
fig_impulse_response = figure('Name', '信道冲激响应分析', 'Position', [150, 150, 1200, 800]);

% 子图1: 冲激响应
subplot(2, 2, 1);
stem(0:length(channel_ir)-1, abs(channel_ir), 'filled', 'Color', colors(1, :));
hold on;
plot(abs(channel_ir), 'Color', colors(2, :), 'LineWidth', 1);
grid on;
xlabel('采样点');
ylabel('幅度');
title(sprintf('信道冲激响应 (%s环境)', channel_type));
legend({'冲激响应', '包络'}, 'Location', 'NorthEast');

% 子图2: 功率延迟分布
subplot(2, 2, 2);
delay_ms = delays * 1000;
powers_db = relative_powers;
stem(delay_ms, powers_db, 'filled', 'Color', colors(1, :));
hold on;
plot(delay_ms, powers_db, 'Color', colors(2, :), 'LineWidth', 2);
grid on;
xlabel('时延 (ms)');
ylabel('相对功率 (dB)');
title('功率延迟分布');

% 子图3: 累积功率分布
subplot(2, 2, 3);
cumulative_power = cumsum(10.^(relative_powers/10));
cumulative_power_norm = cumulative_power / max(cumulative_power);
plot(delay_ms, cumulative_power_norm, 'LineWidth', 2, 'Color', colors(1, :));
hold on;
plot([0, delay_spread*1000], [0.5, 0.5], '--', 'Color', colors(2, :));
plot([delay_spread*1000, delay_spread*1000], [0, 0.5], '--', 'Color', colors(2, :));
grid on;
xlabel('时延 (ms)');
ylabel('累积功率');
title('累积功率分布');
text(delay_spread*1000*1.2, 0.5, sprintf('RMS时延扩展\n%.2f μs', delay_spread*1e6));

% 子图4: 不同环境对比
subplot(2, 2, 4);
channel_types = {'urban', 'suburban', 'indoor', 'rural'};
colors_env = colors(1:length(channel_types), :);

% 串行计算不同环境的参数（避免嵌套并行）
for env_idx = 1:length(channel_types)
    % 根据环境类型设置参数
    switch lower(channel_types{env_idx})
        case 'urban'
            max_delay_env = 2e-3;
            delay_spread_env = 0.5e-3;
            relative_powers_env = linspace(0, -20, num_paths);
        case 'suburban'
            max_delay_env = 1e-3;
            delay_spread_env = 0.2e-3;
            relative_powers_env = linspace(0, -15, num_paths);
        case 'indoor'
            max_delay_env = 0.5e-3;
            delay_spread_env = 0.05e-3;
            relative_powers_env = linspace(0, -12, num_paths);
        case 'rural'
            max_delay_env = 0.2e-3;
            delay_spread_env = 0.02e-3;
            relative_powers_env = linspace(0, -8, num_paths);
    end
    
    % 生成随机时延
    delays_env = sort(rand(1, num_paths) * max_delay_env);
    delays_env(1) = 0;
    
    % 绘制结果
    stem(delays_env*1000, relative_powers_env, 'filled', 'Color', colors_env(env_idx, :));
    hold on;
end
grid on;
xlabel('时延 (ms)');
ylabel('相对功率 (dB)');
title('不同环境对比');
legend(channel_types, 'Location', 'NorthEast');

fprintf('信道冲激响应分析完成！\n');

end